A simple and efficient approach to the optimization of correlated wave functions 
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We present a simple and efficient method to optimize within energy minimization the deter- 
minantal component of the many-body wave functions commonly used in quantum Monte Carlo 
calculations. The approach obtains the optimal wave function as an approximate perturbative solu- 
tion of an effective Hamiltonian iteratively constructed via Monte Carlo sampling. The effectiveness 
of the method as well as its ability to substantially improve the accuracy of quantum Monte Carlo 
calculations is demonstrated by optimizing a large number of parameters for the ground state of 
acetone and the difficult case of the l^Biu state of hexatriene. 



Over the last decade, quantum Monte Carlo (QMC) 
methods have been employed to accurately compute elec- 
tronic properties of large molecular and solid systems 
where conventional quantum chemistry approaches are 
extremely difficult to apply . A crucial step in both the 
variational (VMC) and the diffusion Monte Carlo (DMC) 
approach is the construction of the trial wave function ^'t 
which is usually chosen of the Jastrow-Slater form, that 
is vE't = where $ is a small expansion in Slater de- 
terminants and S the positive Jastrow correlation factor. 

Although considerable progress has been made (princi- 
pally using the variance minimization approach 0) in the 
construction of optimal Jastrow factors, relatively little 
attention has been given to the optimization of the de- 
terminantal part of the wave function. Methods such as 
Hartree-Fock (HF) or a small scale configuration interac- 
tion (CI) are used as a practical way of constructing the 
determinantal component, which is generally not reopti- 
mized when the Jastrow factor is added. However, the 
determinantal part of the wave function solely determines 
the DMC energy Q and often needs to be reoptimized 
to obtain accurate results 0. A practical and simple 
approach to calculate the optimal determinantal compo- 
nent is therefore particularly important to a wide and 
successful application of QMC methods. 

In recent years, several methods to optimize the wave 
function through energy minimization have been pro- 
posed 0, 0, i, Q, S i, IB III III III 111 III 113 A di- 
rect approach to energy minimization entails to compute 
the gradient and the Hessian of the energy with respect 
to the desired parameters. The use of an estimate of 
the Hessian characterized by reduced statistical fluctua- 
tions yields a simple and robust optimization algo- 
rithm for the Jastrow parameters. However, the Hessian 
with respect to the orbital parameters in the determinant 
is affected by higher statistical noise so that devising a 
stable energy-minimization scheme is more difficult hOl: 
the resulting approach is for instance less stable that the 
simple stochastic reconfiguration (SR) method 0, 0| 
and, during the optimization, may have to reduce to the 
inefficient SR to retain stability. To date, the most suc- 
cessful method remains the e nerg y fluctuation potential 
(EFF) method 0, E El E llol which determines the 



optimal determinantal component as the solution of an 
effective Hamiltonian iteratively constructed via Monte 
Carlo sampling. Th e ap proach has been used to opti- 
mize the orbitals 0, llJ , the linear coefficients in front 
of the determinants |j, [l^, and has been extended to 
excited states Q. The method is very stable and more 
efficient than the SR approach but quite complex : the 
construction of the starting effective Hamiltonian as well 
as its update are computationally very demanding steps, 
especially for large systems. 

In this Letter, we propose a simple and efficient op- 
timization method for the determinantal component of 
the wave function. The approach constructs the optimal 
wave function via a perturbative scheme based on the 
EFF method, and only requires easily accessible quanti- 
ties from the quantum chemical calculation used to set 
up the starting wave function. The performance of the 
method is demonstrated on acetone and on the difficult 
case of the l^Bi^ state of hexatriene. First, we will briefly 
review the EFF method. 

The energy fluctuation potential method. Let us as- 
sume that the trial wave function ^'t depends on a set 
of parameters {au}- The derivatives of the energy with 
respect to the parameters can be written as: 



dE 
dak 



= 2{{E^ - E){Ok - Ok)) , 



(1) 



where (•) denotes the average with respect to the square 
of the trial wave function 1^*1 which can be computed 

We dcflned E ^ {Ei^) and 



by Monte Carlo sampling. 
Ok = (Ok) where 



E,^ = 



and 



Ok. = 



1 d-^q 

dak 



(2) 



For the optimal parameters, the derivatives of the energy 
(Eq. ^ are zero, and the fluctuations of the local energy 
i?L and of the functions Ok are uncorrelated. This means 
that the energy is stationary if the remaining fluctuations 
of the local energy cannot be further reduced by adding 
some combination of the functions Ok- Hence, the min- 
imization of the energy can be reformulated as a least- 
squares fit of the fluctuations of the local energy with the 
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functions Ok'- 

x'^m-vo-Y^VkOkr), (3) 

fc>0 

A minimization of with respect to the parameters 14 
leads to the following set of linear equations for m > 0: 

{AEAOrn) = J2 Vl-(AOfcAO,„) (4) 

fc>0 

where we have eliminated Vq from the other equations, 
and AE = Ei, — E and AOm = Om — Om ■ Since the left- 
hand side of these equations corresponds to the deriva- 
tives of the energy (Eq. , the fitting parameters 14 are 
zero if and only if all the derivatives of the energy are 
zero. In general, the parameters 14 which solve these 
linear equations will not be equal to zero and can be 
used to improve the current trial wave function ^^t. 

We focus here on the EFP procedure to optimize the 
determinantal part $ of the trial wave function 

i 

where a spin-adapted configuration state function (CSF) 
Ci is a linear combination of Slater determinants Di of 
single-particle orbitals. Let us assume that the starting 
$ is the lowest solution $q°'' of the CI Hamiltonian 7^^°-': 

7^("'=5:£;f)|<i>|"))($f)|. (6) 

i 

where the states $1"'' span the same space as the CSFs. 
To obtain the optimal coefficients Ci or equivalently the 
best solution in the basis of the eigenstates of Ti.'^^\ we 
consider the variations of $ with respect to the eigen- 
states other than $o"'': 

$ = $(°) ^ = ^4<fi°\ (7) 

fc>0 

so that Ok = ^l^"^ The quantities appearing in the 
linear equations (Eq. [Tj) are sampled from 'fx = J'^o"'' 
and parameters V^^"^ are used to construct a new Hamil- 
tonian Ti.'-^^ as 

^(i)=H(o)+^^(o)(|$(o)^^^(o)|_^^_^^) . (8) 

fc>0 

This Hamiltonian is diagonalized, yielding a new set of 
states, and the procedure is iterated until convergence. 

If we also want to optimize the orbitals in the Slater 
determinants, we can linearize the problem using a so- 
called super-CI expansion to treat the CI and orbital 
variations on the same footing: the CSFs occupied in 
the wave function (Eq. [SJ are augmented by all possible 



single excitations from the occupied to a set of external 
orbitals, and the occupied orbitals can be improved by 
using the natural orbitals of the CI wave function in this 
augmented space The EFP scheme to optimize the 
determinantal component via a super-CI approach 4] is 
quite complex. Setting up the starting super-CI Hamilto- 
nian is not trivial, and the number of singly-excited CSFs 
increases with the number of CSFs in the reference wave 
function and with the size of the basis (i.e. with the num- 
ber of available external orbitals) so that the approach 
quickly becomes computationally too demanding. 

A simplified EFP approach. The basic idea of the 
scheme proposed here is to avoid the explicit construc- 
tion of the EFP Hamiltonian by solving the problem per- 
turbatively. For the example given above, we obtain an 
improved wave function (Eq. [T)) not by diagonalizing the 
Hamiltonian (Eq. IHJ but approximately, to first order in 
the perturbation given by the corrections V^^^ : 

T/(0) 

"i-O ^ ^0 (0) p(0)^fe ■ 

k>0 ~ ^0 

It is simple to show that, in the absence of the Jas- 
trow factor, this is indeed the perturbative solution for 
n-n^°^ since V^°^ = This pertur- 

bative EFP scheme can also be viewed as a generaliza- 
tion of the SR method Tl', T^l which yields an improved 
wave function by applying the operator A — 7i to the 
current state and projecting the result onto the space 
defined by the parameterization. The SR wave func- 
tion is constructed as in Eq. ^ but all the changes are 
scaled by the same energy denominator while the EFP 
approach scales each correction appropriately, achieving 
significantly faster convergence |^. 

Estimate of the energy denominators. The perturba- 
tive solution of the EFP problem requires the knowledge 
of the energy denominators appearing in Eq. ^ which 
can be easily estimated with the use of quantities readily 
available from the quantum chemical calculation used to 
set up the starting wave function. We will first consider 
the optimization of the orbital parameters. 

The simplest case is given by a trial wave function 
constructed from a single closed-shell determinant as 
^I^T = JDq, where Dq is a determinant of Hartree-Fock 
or density functional Kohn-Sham orbitals which we want 
to reoptimize in the presence of the Jastrow factor. We 
assume here that we have chosen HE orbitals but all the 
following considerations equally apply to the case of den- 
sity functional orbitals. The set of occupied and virtual 
HE orbitals and the atomic basis on which they are ex- 
panded span the same space so that we can express the 
variations of one orbital with respect to the expansion 
coefficients in terms of the variations with respect to the 
other orbitals. To optimize the M occupied orbitals Lpi in 
Dq, we only need to consider the variations with respect 
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to the N — M virtual orbitals as 



N 

j=M+l 



(10) 



where we only mix orbitals of the same symmetry. The 
corresponding first-order change in the wave function is 



M N 

E E 

i=l j=M+l 



(11) 



where Cq"*"* — Dg"^"'" 



D. 







is the CSF of the two 



determinants obtained by substituting the orbital ipi with 
orbital Lpj for the up and down spin, respectively. It is 
simple to show that 0'~*^ = Cq^"'/Z?o- 

Following the EFP procedure, we sample from the 
quantities appearing in the linear equations (Eq. and 
obtain the corrections V^^^ corresponding to the func- 
tions O'^^ . We can then obtain an improved determi- 
nantal component in analogy to Eq. as 



M 



N 



*'-^o-E E 



=1 i=A/+i 



(12) 



where we assign an energy scale AE^^^ to the varia- 
tion corresponding to the single excitation Cq^-*. We 
now make the key observation that, to first order, this is 
equivalent to the much simpler step of constructing an 
improved set of occupied single-particle orbitals as 



<y5j = '^t 



N 

E xe 



(13) 



which can be used in a new single-determinant wave func- 
tion — JD'q. We then proceed iteratively by sam- 
pling a new set of V^^^ from the wave function ^'r^ and 
updating the orbitals as in Eq. (|13|l until convergence. 

To estimate the energy denominators /S.E'^^^ , we note 
that the starting HF orbitals are canonical orbitals ob- 
tained by diagonalizing the Fock matrix so that the cor- 
responding eigenvalues have a physical interpretation via 
the Koopman's theorem. Minus the eigenvalues of the oc- 
cupied and the virtual orbitals are an approximation to 
the ionization potentials (IP) and the electron affinities 
(EA), respectively. The energy associated to the pro- 
motion of one electron from orbital (pi to orbital ipj can 
therefore be estimated as the difference IP— EA, namely 



AE' 



HF 



(14) 



where and are the eigenvalues corresponding to or- 
bitals Lpi and Lpj. We do not find it necessary at succes- 
sive iterations to update the parameters AE'^^^ which 
are therefore kept fixed at the starting HF values. 

In Fig. ^ we compare the convergence of the proposed 
method, denoted by "perturbative EFP" , with the EFP 
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FIG. 1: Convergence of the VMC energy of acetone with the 
perturbative EFP, the EFP and the SR method to optimize 
the 473 orbital parameters. The starting determinantal com- 
ponent is from an unconverged HF calculation. The statistical 
errors are smaller than the symbol size. 



and the SR methods for the orbital optimization of a 
single closed-shell Jastrow-Slater wave function for the 
ground state of acetone (CaHgO). To demonstrate the 
robustness of the method, we construct a poor starting 
determinantal component from a set of unconverged HF 
orbitals obtained by diagonalizing the Fock matrix built 
from the Huckel orbitals. A simple electron-electron and 
electron-nucleus Jastrow factor is employed and the 473 
orbital parameters are optimized using the starting un- 
converged eigenvalues as estimates for the energy denom- 
inators (Eq. I14|l . The SR optimization is performed with 
the critical value of A estimated as in Ref. 0, and con- 
verges to the optimal wave function only after 40 steps. 
The perturbative EFP approach converges in only 2 steps 
as the EFP method, with a computational effort as low as 
the one of the SR method. Interestingly, the perturbative 
method is successful even when using the eigenvalues of 
the starting Huckel orbitals to roughly estimate the en- 
ergies AE^^^ . Finally, if we start from converged HF 
orbitals, the energy gain from the optimization is of 14 
± 1 and 6 ± 1 mHartree in VMC and DMC, respectively. 

The procedure is not very different if the wave function 
is constructed from a spin-unrestricted determinant. The 
variations of an occupied orbital fEa. 110(1 are now with re- 
spect to the virtual orbitals of the same spin, and the ex- 
citation energies are estimated as AE'^'"'^^-°' ^ ej^a — ^i.a 
where a is the spin of the varied orbital. The case 
of a restricted open-shell determinant is instead more 
complicated since the final orbitals and eigenvalues of 
a restricted open-shell Hartree-Fock calculation are not 
uniquely defined. Among the many constructions avail- 
able in the literature, we find that the one by Guest and 
Sanders [l^ yields eigenvalues which, in average, well ap- 
proximate the excitation energies in Eq. ((14|l . Note that 
closed-to-open-shell excitations must also be included in 
Eq. H10|) to allow full orbital variation. 

For a multideterminant wave function from a complete- 
active-space self-consistent-field (CASSCF) calculation, 
a set of canonical orbitals is usually defined as diagonal- 
izing the closed, active and virtual blocks of the so-called 
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State 


Wave function 


EvMC 


Edmc 


AE (eV) 




HF 


-38.684(1) 


-38.7979(7) 


- 




B3LYP 


-38.691(1) 


-38.7997(7) 


- 




optimized 


-38.691(1) 


-38.7992(7) 


- 


l'Bi„ 


CAS(2,2) 


-38.472(1) 


-38.5910(7) 


5.63(3) 




B3LYP 


-38.482(1) 


-38.6030(7) 


5.35(3) 




optimized 


-38.493(1) 


-38.6069(8) 


5.23(3) 


expt. [ 


20] 






5.22 



TABLE I: VMC and DMC energies in Hartree for the l^Ag 
and the l^Biu states of trans-hexatriene. The DMC excita- 
tion energies AE are with respect to the ground state ob- 
tained with the same wave function type. We optimize 604 
and 493 parameters for the ground and excited state, respec- 
tively. The errors on the last figure are given in parenthesis. 



generalized Fock matrix. However, to obtain a good es- 
timate of the closed-shell-to-active and active-to- virtual 
excitation energies, the eigenvalues must be properly ad- 
justed since the energy of an active orbital should be 
different for excitations into or out of it, i.e. closer to 
an electron affinity or to an ionization potential, respec- 
tively. Therefore, for these excitations, we follow Ref. 
to define AE^^^ ~ E^^^ - i^cAsscF as: 



AE' 



+ t; (Pn + 2 - Pu) (15) 



where e is a CASSCF eigenvalue, p is the single-particle 
density matrix, and A is an average difference between 
the electron affinity and the ionization potential of the 
active orbitals, chosen of the order of 0.3-0.5 Hartree. 

Finally, we consider the simultaneous optimization of 
orbitals and CI coefficients. The simplest approach is to 
alternate between an orbital optimization step with the 
perturbative EFP method, and a CI calculation in the 
basis of the CSFs multiplied by the Jastrow factor [||. 
Alternatively, we can use the perturbative EFP method 
for orbital and CI parameters and sample the quantities 
in Eq. Q for both variations. The corrections 14 for the 
CI coefficients are computed in the basis of the CI states 
obtained in a starting CI calculation, and the energy de- 
nominators (Eq. O estimated using the CI energies. 

In Table U we demonstrate the performance of the 
method on the l^Bi„ state of trans-hexatriene (CgHg) 
which represents a challenge for all electronic structure 
approaches: CASPT2 yields a vertical excitation energy 
of only 5.01 eV j21j while the TDDFT energies range 
between 4.42 and 4.64 eV, depending on the functional 
employed ^2^. For the single-determinant closed-shell 
ground state, the VMC energy is improved by about 7 
mHartree when otimizing the HF orbitals while, differ- 
ently from acetone, the DMC energy is not significantly 
affected [i^. The excited state is mainly an excita- 
tion from the highest occupied (LUMO) to the lowest 
unoccupied (HOMO) orbital and well described a two- 
determinant singlet wave function. A CASSCF wave 



function of two electrons in two orbitals, denoted by 
CAS(2,2), gives a DMC vertical excitation energy 0.4 eV 
higher than the experimental value. Employing B3LYP 
orbitals partially improves the result but only the fully 
optimized wave function yields an excitation energy in 
perfect agreement with experiments. Finally, we note 
that both ground and excited states do not have multi- 
configurational character and using an unoptimized mul- 
tideterminant wave function yields significantly worse 
QMC energies: for the excited state, a CASSCF wave 
function of 6 electrons in 6 orbitals, where we keep the 
leading 5 CSFs, yields an energy of -38.438(1) and - 
38.560(2) Hartree in VMC and DMC, respectively Upon 
optimization of both orbitals and CI coefficients, we re- 
cover the best energies of Table 
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